Add geometric multigrid Poisson solver (scalar backend)#200
Open
shaia wants to merge 4 commits into
Open
Conversation
CG iteration counts grow O(N) with grid refinement; geometric multigrid converges in a grid-size-independent number of cycles for O(N) total work, making it the optimal solver for large structured-grid pressure Poisson problems (ROADMAP 1.2). V/W/F(FMG) cycles with Red-Black Gauss-Seidel or weighted-Jacobi smoothers, full-weighting restriction and bilinear/trilinear prolongation, 2D/3D. Grid dims must be 2^k+1 per active dimension. Two BC modes, because the spec's Dirichlet assumption did not match the codebase: MG_BC_NEUMANN (default) mirrors the zero-gradient convention of the other solvers, MG_BC_DIRICHLET holds caller boundary values fixed. The Neumann mirror slaves the boundary vertex to its interior neighbor, so restriction folds the toward-boundary weights to stay the exact adjoint of prolongation - without the fold the convergence factor degrades several-fold. The singular Neumann system is handled by mean-projecting restricted RHS vectors and coarse corrections; the coarsest Neumann grid stops at 5 points per dim (a 3x3 mirror grid has an identically zero operator). Neumann cross-checks use RB-SOR, not CG: CG freezes boundary values during Krylov iterations (BC applied only at start/end) and therefore solves a different discrete system than the per-sweep-BC solvers. Also fix poisson_solve_3d ignoring init status - a failed init (now possible via the 2^k+1 check) no longer caches a broken solver. Verified: V(2,2) factor < 0.15, 9x9-129x129 converge in <= 15 cycles (spread <= 3), FMG reaches discretization accuracy in one pass; full suite 103/103 incl. Ghia/cavity validation.
Contributor
There was a problem hiding this comment.
Pull request overview
This PR adds a new geometric multigrid (GMG) Poisson solver targeting the scalar CPU backend, aiming to deliver grid-size-independent convergence and O(N) work on structured grids, along with operator validation, convergence tests, and documentation updates.
Changes:
- Implement geometric multigrid scalar solver with V/W/F(FMG) cycles, Neumann/Dirichlet BC modes, and 2^k+1 dimension validation.
- Add multigrid transfer operators (restriction/prolongation) plus dedicated operator and convergence/unit tests.
- Extend Poisson solver API/docs/roadmap/changelog and add a new convenience preset + cache-safety fix in
poisson_solve_3d().
Reviewed changes
Copilot reviewed 16 out of 16 changed files in this pull request and generated 2 comments.
Show a summary per file
| File | Description |
|---|---|
| tests/math/test_multigrid_operators.c | New unit tests validating restriction/prolongation operators, adjointness, conservation, and mean handling. |
| tests/math/test_multigrid_convergence.c | New end-to-end solver tests for convergence factors, grid independence, BC modes, and API/caching behavior. |
| ROADMAP.md | Marks GMG as implemented and documents the 2^k+1 grid-size limitation. |
| lib/src/solvers/linear/multigrid_internal.h | Introduces internal MG structs/utilities and declares transfer operators for internal/test use. |
| lib/src/solvers/linear/linear_solver.c | Wires MULTIGRID into factory creation + adds convenience preset caching and init-failure handling. |
| lib/src/solvers/linear/linear_solver_internal.h | Declares the multigrid scalar factory function for internal linkage. |
| lib/src/solvers/linear/cpu/multigrid_transfer.c | Implements full-weighting restriction, (tri)linear prolongation, and interior-mean helpers. |
| lib/src/solvers/linear/cpu/linear_solver_multigrid.c | Implements the scalar multigrid solver lifecycle and cycle logic (V/W/FMG), smoothers, residuals, transfers. |
| lib/include/cfd/solvers/poisson_solver.h | Adds MG enums/params, solver type constant, and convenience enum for MG scalar. |
| lib/CMakeLists.txt | Adds multigrid sources to the scalar build. |
| examples/poisson_solver_tuning.c | Updates example to demonstrate “unavailable backend” behavior for multigrid on GPU. |
| docs/reference/solvers.md | Documents GMG algorithm, constraints, parameters, and adds to performance comparison. |
| docs/reference/api-reference.md | Updates API reference to include multigrid method and parameters. |
| docs/guides/examples.md | Updates example guide text to reflect multigrid backend availability constraints. |
| CMakeLists.txt | Adds new test executables and registers them with CTest. |
| CHANGELOG.md | Adds release notes for new multigrid solver and the init-status caching fix. |
💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.
The mg_* operators are already part of CFD::Library in static builds (the default here), so compiling multigrid_transfer.c into the test binary risks duplicate definitions; only shared builds, where the symbols are hidden, need the direct compilation. Addresses PR #200 review feedback.
Callers that only see NULL (e.g. the error-handling demo reading cfd_get_last_status) previously got a stale status. Report CFD_ERROR_UNSUPPORTED for unavailable backends and CFD_ERROR_INVALID for unknown methods, matching the factory convention in cfd_solver_create. The non-multigrid branches had the same pre-existing gap and are fixed uniformly. Addresses PR #200 review feedback.
The operators are not declared static — they have external linkage and are merely not exported from the shared library. Saying "internal linkage" misstated how the unit test is able to build against them. Addresses PR #200 review feedback.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
CG iteration counts grow O(N) with grid refinement; geometric multigrid converges in a grid-size-independent number of cycles for O(N) total work, making it the optimal solver for large structured-grid pressure Poisson problems (ROADMAP 1.2).
V/W/F(FMG) cycles with Red-Black Gauss-Seidel or weighted-Jacobi smoothers, full-weighting restriction and bilinear/trilinear prolongation, 2D/3D. Grid dims must be 2^k+1 per active dimension.
Two BC modes, because the spec's Dirichlet assumption did not match the codebase: MG_BC_NEUMANN (default) mirrors the zero-gradient convention of the other solvers, MG_BC_DIRICHLET holds caller boundary values fixed. The Neumann mirror slaves the boundary vertex to its interior neighbor, so restriction folds the toward-boundary weights to stay the exact adjoint of prolongation - without the fold the convergence factor degrades several-fold. The singular Neumann system is handled by mean-projecting restricted RHS vectors and coarse corrections; the coarsest Neumann grid stops at 5 points per dim (a 3x3 mirror grid has an identically zero operator).
Neumann cross-checks use RB-SOR, not CG: CG freezes boundary values during Krylov iterations (BC applied only at start/end) and therefore solves a different discrete system than the per-sweep-BC solvers.
Also fix poisson_solve_3d ignoring init status - a failed init (now possible via the 2^k+1 check) no longer caches a broken solver.
Verified: V(2,2) factor < 0.15, 9x9-129x129 converge in <= 15 cycles (spread <= 3), FMG reaches discretization accuracy in one pass; full suite 103/103 incl. Ghia/cavity validation.